Thermodynamics of the 3D Hubbard model on approach to the Neel transition 
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We study the thermodynamic properties of the 3D Hubbard model for temperatures down to the 
Neel temperature using cluster dynamical mean-field theory. In particular we calculate the energy, 
entropy, density, double occupancy and nearest-neighbor spin correlations as a function of chemical 
potential, temperature and repulsion strength. To make contact with cold-gas experiments, we also 
compute properties of the system subject to an external trap in the local density approximation. 
We find that an entropy per particle S/N ~ 0.65(6) at U/t = 8 is sufficient to achieve a Neel state 
in the center of the trap, substantially higher than the entropy required in a homogeneous system. 
Precursors to antiferromagnetism can clearly be observed in nearest-neighbor spin correlators. 
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The Hubbard model [lj remains one of the corner- 
stone models in condensed matter physics, capturing 
the essence of strongly correlated electron physics rele- 
vant to high-temperature superconductors [2| and corre- 
lation driven insulators 0]. While qualitative features 
of the phase diagram are known from analytical ap- 
proximations, controlled quantitative studies in the low- 
temperature regimes relevant for applications are not 
readily tractable with tools presently available. A recent 
program that aims to implement the Hubbard model in a 
cold gases experiment 4 1 has led to experimental signs of 
the Mott insulator M. Model ing by dynamical mean 
field theory (DMFT) 00] and high-temperature series 
expansions [8[ resulted in temperature and entropy es- 
timates and an error budget [9J. A major experimental 
achievement will be the detection of the antiferromag- 
netic phase, for which the slow and ill-understood equili- 
bration rates, the limited number of detection methods, 
and inherent cooling problems will have to be overcome. 

Experimental progress has also sparked interest in sim- 
ulations of the 3D Hubbard model where new algorithms 
have been developed to treat the Hubbard model, such 
as the real-space DMFT [HI El or diagrammatic Monte 
Carlo [lH studies. Similar to the case of bosons, where 
synergy between experiment and simulation has led to 
quantitative understanding of experiments fl3| , accurate 
results for the thermodynamics of the 3D Hubbard model 
will also be crucial for validation, calibration and ther- 
mometry of fermionic experiments. A crucial role is 
played by the entropy, since these experiments form iso- 
lated systems where ideally the parameters are changed 
adiabatically (and not isothermally) . 

In this Letter we provide the full thermodynamical 
equation of state of the Hubbard model - in particu- 
lar the entropy, energy, density, double occupancy and 



spin correlations - for interactions U up to the band- 
width 12t on approach to the Neel temperature Tjv by 
performing controlled large-scale cluster dynamical mean 
field calculations and extrapolations to the infinite sys- 
tem size limit, as well as dctcrminantal diagrammatic 
Monte Carlo (DMC) simulations at half filling. We then 
use this information to calculate the entropy per particle 
required for experiments on ultracold atomic gases in op- 
tical lattices to reach a Neel state in the trap center. We 
finally show that the nearest-neighbor spin correlation 
function contains clear precursors for antiferromagnetism 
that may already be detectable in current generation ex- 
periments and that are useful for thermometry (more so 
than measurements of the double occupancy) close to T/v. 
The Hubbard model is defined by its Hamiltonian 



H 



U^fn 



fafa 



(1) 



where c| CT creates a fermion with spin component a =t, 4- 
on site i, hi a = c\ a c ia , (...) denotes summation over 
neighboring lattice sites, t is the hopping amplitude, U 
the on-site repulsion, and fa = fj, — V(fi) with \x the 
chemical potential and ~V(fi) the confining potential at 
the location of the i-th lattice site. We will set V(r) = 
in all calculations and consider realistic traps later on. 

Our numerical approach is a cluster generalization of 
dynamical mean field theory [14J. In cluster DMFT the 
self energy is approximated by iV c momentum-dependent 
basis functions ^(fc): cj) ps J2k < t } K{k)Yj K (uj). The 
exact problem is recovered for N c —> oo. Within the 
dynamical cluster approximation (DCA)[l5| used here, 
4>K{k) are piecewise constant over momentum patches. 
The DMFT method (3| is the N c = 1 cluster approxi- 
mation where £ = S(w) and no momentum dependence 
is retained. 
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FIG. 1: (color online) Entropy per lattice site s of the Hub- 
bard model as a function of temperature T/t, for U/t = 8, at 
half filling. Shown with dashed vertical lines is Tjv from Ref. 
[2l| (black) and with dotted lines according to DMC (blue). 
Shown with dashed horizontal lines is the entropy per lattice 
site s at Tjv [2l|] (black), with dotted lines according to DMC 
(blue); also log(2) is shown as a full horizontal line (green). 
The inset shows the energy E per lattice site. 



Solving the DMFT and DCA equations requires the 
solution of a quantum impurity model. Continuous- 
time quantum impurity solvers [17H19I ] , in particular the 
continuous-time auxiliary field (CT-AUX) method used 
here 18( , combined with advanced numerical techniques 
[20| , have made it feasible to solve such models efficiently 
and numerically exactly on large clusters, thereby pro- 
viding a good starting point for an extrapolation of finite 



size clusters to the infinite system [12], |21|, |22| . We have 
performed extensive DCA calculations on bipartite clus- 
ters with N c = 18,26,36,48,56, and 64. In order to 
achieve an optimal scaling behavior we exclusively use 
the clusters determined in Ref. [2l[ following the criteria 
proposed by 23j . As the DCA exhibits a 1/L 2 finite-size 

1 /3 I — 

scaling in the linear cluster size we ex- 

trapolate our cluster results linearly in N c . Our error 
bars include extrapolation uncertainties. Despite a sign 
problem away from half-filling, temperatures T/t > 0.4, 
on the order of the Neel temperature, are reliably acces- 
sible for all but the largest interaction strength U = 12t 
where we have been restricted to T/t > 0.5. 

The potential energy, double-occupancy, and nearest- 
neighbor spin-spin correlation have been measured di- 
rectly. The kinetic energy E^ n = J2 n % e(k)G(k, iuj n ) has 

been calculated by summing e(k), the bare-dispersion of 
the simple cubic lattice, and the single-particle Green 
function G(k,iui n ) over all momenta k and Matsubara 
frequencies iuj n . The entropy S has subsequently been 




FIG. 2: (color online) Double occupancy of the Hubbard 
model as a function of temperature T, at half filling. Ex- 
trapolated DCA results are shown as solid lines and DMFT 
values as dashed lines. The vertical lines are the same as in 
Fig.ffl 



calculated by numeric integration 



S(T) = S(T U ) 



E(T U ) + E(T) 



T, 



T 



Tu dT >^n 

, a T a : 



(2) 



up to a T u /t ps 10 , where the entropy S(T U ) is accurately 
given by a high-temperature series expansion. Tables of 
the complete results containing finite cluster and extrap- 
olated values at and away from half filling for the entropy, 
energy, density, double occupancy, and spin correlations 
are given in the supplementary material [251 ] . 

Results at half filling - We start our analysis at half 
filling and focus on U/t = 8 where we can compare with 
results from DMC simulations (2(| (and the more accu- 
rate [27]]) ■ We see in Fig. Q] that the entropy calculated 
using DCA and DMC coincides within error bars at all 
temperatures. Agreement with a 10th order high tem- 
perature series expansion Q is found down to T/t ps 1.6. 
At that temperature also single site DMFT starts to de- 
viate because that method misses short-range antiferro- 
magnetic correlations. The Neel temperature was found 
to be T N /t ps 0.36(2) in Ref. ED. Our DMC calculations 
find it more accurately at T^/t = 0.333(7). Using DCA 
calculations, the critical entropy is s ps 0.46(2) for Tm 
according to Ref. [HJ, and s := S/N c ps 0.41(3) with T N 
according to the DMC. In the rest of the paper we will 
only use the Tjv as determined by DMC but with en- 
tropies calculated by DCA (since away from half filling 
only entropies calculated by DCA are available). 

The double occupancy, which hasplayed a crucial role 
in optical lattice experiments [H, m, B 0] , is shown in 
Fig. [2] as a function of temperature at half filling for dif- 
ferent values of U/t. While for small U/t a remarkable 
increase is seen on approach to Tjv, only a plateau re- 
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FIG. 3: (color online) Nearest-neighbor spin-spin correlation 
of the Hubbard model as a function of temperature T, at half 
filling. The inset shows the density dependency for U/t = 8 
and selected temperatures. Vertical lines: compare to Fig. \T\ 
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FIG. 4: (color online) Entropy per lattice site s and the en- 
tropy per particle S/N (inset) of the Hubbard model at the 
temperature T = 0.35t ~ Tjv, as a function of density n, for 
U/t — 8. DMFT values and extrapolated DCA values are 
shown. The errors are dominated by extrapolation errors. 



mains at moderate values of U/t. This is in contrast 
to the DMFT predictions, but similar to lattice QMC 
results in two dimensions (30j . For larger interactions 
(U/t > 12), the double occupancy rises above that of 
a single site paramagnet, consistent with DMFT results 
for the anti- ferromagnetic phase below T/v ll| . The neg- 



ative slope of D(T), discussed in the context of single 
site DMFT [28|, persists for a wide range of parameters. 
Sharp features just above Tjv, as detected in single site 
(momentum independent) studies ll|, are not observed 
for the interaction values and temperature ranges studied 
here. Hence the proposal that the double occupancy is a 
good candidate for thermometry is not substantiated by 
more accurate momentum-dependent calculations. 

The spin-spin correlation function plotted in Fig. [3] as 
a function of temperature for various U and as a function 
of filling for T/t = 0.3,0.8, and 1.6 at U/t = 8, is only 
accessible in methods that include non-local correlations, 
but may be accessible experimentally (3l| . It has a steep 
slope on approach to the Neel temperature, which makes 
it an ideal quantity for thermometry. This corresponds 
to the intuitive picture that charge degrees of freedom 
are already essentially frozen out around T/v while the 
spin degrees of freedom start to order there. 

Results away from half filling - Fig.[4]shows the entropy 
per lattice site for U/t = 8. The inset demonstrates that 
entropy per particle number N increases strongly at lower 
densities. While single site DMFT remains accurate for 
densities lower than n < 0.6 due to the weak momen- 
tum dependence of the self-energy in this regime [HI, 
the DCA results are important for larger densities. Sim- 
ilarly, near half filling DMFT overestimates the double 
occupancy by 10%, while deviations are less pronounced 
at lower densities. This observation persists for all in- 



teractions and temperatures investigated. The flattening 
of the double occupancy for 8 < U/t < 12 (cf. Fig. [2]) 
is also seen away from half filling, and leads to virtually 
unchanged profiles over the trap in an optical lattice sys- 
tem. On the other hand, the spin-spin correlation func- 
tion away from half filling (see the inset in Fig. [3]) changes 
most rapidly near half filling when approaching the Neel 
temperature since it couples strongly to the developing 
(short-range) spin correlations. 

Entropy in the optical lattice system — We now turn to 
the experimentally relevant case of an optical lattice in 
a harmonic trap, which is a closed system where entropy 
is conserved and temperature changes when adiabati- 
cally changing the parameters of the Hamiltonian. Wc 
choose parameters close to current experiments: V(r) = 
0.004(|r|/a) 2 i with lattice spacing a, and we consider the 
case of half filling in the trap center: fi = U/2. We treat 
the harmonic confinement in a local density approxima- 
tion(LDA): for every site we perform a DCA simulation 
for a homogeneous system and average the results over 
the trap. LDA was found to be a good approximation for 
the Bosc-Hubbard model for wide traps, except in close 
proximity to the critical point 33-3{| of the U(l) phase 
transition because of the diverging correlation length and 
in our setup errors due to the LDA can be neglected com- 
pared to our other systematic errors. 

Due to the large volume fractions, the wings of the 
gas may capture more entropy than the center of the 
trap, even though the entropy per site is comparable to 
the one in the center (see Fig. [5]). In fact, the entropy 
of the whole density range 0.1 < n < 0.9 is large, and 
this opens the possibility to observe anti-ferromagnetic 
order in the trap center at an average entropy per par- 
ticle over the trap which is about 50% larger than what 
could be expected from a homogeneous study. Optimal 
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FIG. 5: (color online) Entropy profiles (entropy per lattice 
site) plotted over the trap in the LDA approximation for dif- 
ferent temperatures with an interaction strength U/t = 8. 
Error bars are shown every 5 lattice spacings, but are gener- 
ally smaller than the symbol size. 

parameters are around U/t — 8 when Tjv/i = 0.333(7) 
according to DMC, corresponding to S/N = 0.65(6) in 
the trap, while S/N = 0.41(3) would be expected for 
a homogeneous system. As seen in Fig. [6j all U in the 
range 8 < U/t < 12 lead to similar conclusions. We have 
verified that changing the trap by a factor of 4 does not 
alter these conclusions. 

Conclusions - We have provided the full thermody- 
namics of the 3d Hubbard model using the DCA for- 
malism for U/t < 12 and temperatures above the Neel 
temperature. Comparing to single site DMFT results we 
found that the latter already fail at remarkably high tem- 
peratures (T/t « 1.5 for U/t = 8 at the 1% level) near 
half filling. While the entropy per particle at the Neel 
temperature T^/t = 0.333(7) (determined with DMC) is 
S/N = 0.41(3) for U/t = 8 in a homogeneously half filled 
system, we find that the Neel transition in a trap can al- 
ready be reached at S/N — 0.65(6) in a realistically sized 
harmonic trap (taking Tjy according to Ref. [2lJ leads to 
S/N = 0.69). 

We have also investigated the double occupancy and 
the nearest-neighbor spin-spin correlation function as 
quantities that are experimentally measurable and which 
were suggested to show precursors of antiferromagnetism. 
It turns out that the double occupancy is more or less flat 
as a function of temperature, while the spin correlations 
show a strong temperature dependence around the Neel 
temperature. This suggests that the spin correlations, 
not the double occupancy, are best suited to observe pre- 
cursors of antiferromagnetism and measure the temper- 
ature. Our numerical data can be used to calibrate such 
a spin-correlation thermometer. 
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FIG. 6: (color online) Entropy per particle averaged over the 
trap as a function of temperature relative to the Neel tem- 
perature for different U. The Neel temperature is reached 
at T/t — 0.333(7) for U/t = 8, when the average entropy is 
S/N = 0.65(6). Errors (not shown) in S(Tjv) are estimated 
to be in the 10% range, with the largest contribution caused 
by the uncertainty in Tn- 
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